##-------------------------------------------------------------------------------##
# Replication code for "Priming political trust: Evidence from an experiment"
# 
# author: kyle peyton (kyle.peyton@yale.edu)
##-------------------------------------------------------------------------------##

require(ggplot2)
require(scales)
require(ggthemes)
require(gridExtra)

load("probity.RData")

#### ----------------------------- Table 1 ------------------------------------ ####
tapply(probity$confgov, probity$trusttreat, function (x) mean(x,na.rm=T))
tapply(probity$confgov, probity$trusttreat, function (x) sd(x,na.rm=T))

tapply(probity$fedtrust, probity$trusttreat, function (x) mean(x,na.rm=T))
tapply(probity$fedtrust, probity$trusttreat, function (x) sd(x,na.rm=T))

tapply(probity$statetrust, probity$trusttreat, function (x) mean(x,na.rm=T))
tapply(probity$statetrust, probity$trusttreat, function (x) sd(x,na.rm=T))

tapply(probity$healthtrust, probity$trusttreat, function (x) mean(x,na.rm=T))
tapply(probity$healthtrust, probity$trusttreat, function (x) sd(x,na.rm=T))

tapply(probity$unitrust, probity$trusttreat, function (x) mean(x,na.rm=T))
tapply(probity$unitrust, probity$trusttreat, function (x) sd(x,na.rm=T))

#### ----------------------------- T-tests ------------------------------------ ####
Z <- probity$trusttreat
t.test(probity$poltrust ~ Z)
t.test(probity$selfserv ~ Z)
t.test(probity$bigint ~ Z)
t.test(probity$corrupt ~ Z)

t.test(probity$fedtrust ~Z)
t.test(probity$statetrust ~Z)
t.test(probity$healthtrust ~Z)
t.test(probity$unitrust ~Z)

#### --------------------- Analysis and Figure 2 ----------------------------- ####
probity$full <- rep(1,nrow(probity))

dvs <- c("confgov","fedtrust","statetrust","unitrust","healthtrust") ## define dvs of interest
conds <- c("full")
trust.estate <- matrix(NA, nrow= length(dvs), ncol = length(conds))
trust.atese <- matrix(NA, nrow= length(dvs), ncol = length(conds))

## fix to handle NAs
fix <- rep(NA,ncol(probity))
probity <- rbind(probity,fix)

for (i in 1: length(dvs)){
    for (j in 1:length(conds)) {
        
        Y <- probity[-which(is.na(probity[,dvs[i]]) | is.na(probity[,conds[j]])),][,dvs[i]]    
        Z <- probity[-which(is.na(probity[,dvs[i]]) | is.na(probity[,conds[j]])),]$trusttreat
        X <- probity[-which(is.na(probity[,dvs[i]]) | is.na(probity[,conds[j]])),][,conds[j]]  
        
        ##standardize the outcome variable for interpretation and sign
        Y <- as.numeric((Y-mean(Y))/sd(Y))*-1
        
        trust.estate[i,j] <- mean(Y[Z==1 & X == 1]) - mean(Y[Z==0 & X == 1])
        trust.atese[i,j] <- sqrt(var(Y[Z==1 & X == 1])/length(Y[Z==1 & X == 1]) + var(Y[Z==0 & X == 1])/length(Y[Z==1 & X == 1]))
        
    }
}

## assign row and column names
rownames(trust.estate) <- dvs
colnames(trust.estate) <- conds

rownames(trust.atese) <- dvs
colnames(trust.atese) <- conds

poltrust <- data.frame(trust.estate[1,], trust.atese[1,], "Political Trust")
fedtrust <- data.frame(trust.estate[2,], trust.atese[2,], "Federal Trust")
statetrust <- data.frame(trust.estate[3,], trust.atese[3,], "State Trust")
unitrust <- data.frame(trust.estate[4,], trust.atese[4,], "University Trust")
healthtrust <- data.frame(trust.estate[5,], trust.atese[5,], "Health Trust")

colnames(poltrust) <- c("te","se","group")
colnames(fedtrust) <- c("te","se","group")
colnames(statetrust) <- c("te","se","group")
colnames(unitrust) <- c("te","se","group")
colnames(healthtrust) <- c("te","se","group")

probitysub <- rbind(poltrust,fedtrust,statetrust,healthtrust,unitrust)

probityplot <- ggplot(probitysub, aes(x=group, y = te, fill = group)) + 
    geom_point(size = 5, shape = 21, fill = "black", colour = "black",lwd=3) +  
    xlab(NULL) + ylab("Standardized Treatment Effect") + 
    geom_hline(yintercept=0, colour = "black", size = 0.5) + 
    geom_errorbar(aes(ymin=te-(1.96*se), ymax=te+(1.96*se)),
                  size=1,
                  width=.5,
                  position=position_dodge(.9)) + 
    ggtitle("") +
    scale_y_continuous(breaks=seq(-1.25,1.25,.25), limits = c(-1.25,1.25)) +
    theme_bw() +
    theme(axis.ticks.x = element_blank(), 
          legend.title=element_blank(),
          legend.position = c(.12,0.15),
          panel.grid.minor.x=element_blank(), 
          panel.grid.major.x=element_blank(),
          panel.grid.minor.y=element_blank(),
          legend.text = element_text(color= 'black', size = 8),
          legend.key = element_rect(fill = 'white', size = 0.5, linetype='solid'),
          legend.key.size = unit(1, 'lines'),
          legend.background = element_rect(fill="transparent"),
          text = element_text(family = "Palatino", size = 20, face = "bold")) + 
    geom_vline(xintercept = 1.5:4.5, colour="grey", linetype = "longdash") 
ggsave(probityplot, file= "inttrustplot.pdf",height = 10, width=17)




#### --------------------- Balance/Analysis for Appendix ----------------------------- ####
require('exactRankTests')

## exact test for pre-treatment covariates 
wilcox.exact(probity$gender~Z,exact = T)
wilcox.exact(probity$coalition~Z,exact = T)
wilcox.exact(probity$green~Z,exact = T)
wilcox.exact(probity$labor~Z,exact = T)
wilcox.exact(probity$ideology~Z,exact = T)
wilcox.exact(probity$age~Z,exact = T)


## balance table
pretreat <- c("gender","ideology","age")
balance <- matrix(NA,ncol=2,nrow=3*length(pretreat))
j=1;k=3
for(i in 1:length(pretreat)){
    balance[j:k,1] <- matrix(data = round(c(mean(probity[probity$trusttreat==1,][[pretreat[i]]],na.rm=T),
                                            sd(probity[probity$trusttreat==1,][[pretreat[i]]],na.rm=T),
                                            sum(is.na(probity[probity$trusttreat==1,][[pretreat[i]]]))),2),nrow=3,ncol=1)
    balance[j:k,2] <- matrix(data = round(c(mean(probity[probity$trusttreat==0,][[pretreat[i]]],na.rm=T),
                                            sd(probity[probity$trusttreat==0,][[pretreat[i]]],na.rm=T),
                                            sum(is.na(probity[probity$trusttreat==0,][[pretreat[i]]]))),2),nrow=3,ncol=1)
j=j+3;k=k+3
}

table(probity$partyid,probity$trusttreat)






